import numpy as np
def zsgrad(fun,gfun,Hfun,x0,k,eps,maxk=100):
    gk= gfun(x0)
    if np.linalg.norm(-gk) < eps or k >= maxk:
        return x0, fun(x0), k
    Hk = Hfun(x0)
    dk = -np.linalg.inv(gk.reshape(len(gk),-1) * Hk * gk) * (gk.reshape(len(gk),-1) * gk)
    return zsgrad(fun,gfun,Hfun,x0+np.dot(dk,gk),k+1,eps,maxk)

result = zsgrad(fun = lambda x:100*(x[1]-x[0])**2 + (1-x[0])**2,
                        gfun = lambda x : np.array([-400*x[0]*(x[1]-x[0]**2)-2*(1-x[0]),200*(x[1]-x[0]**2)],'float32'),
                        Hfun = lambda x : np.array([[1200*(x[0]**2)-400*x[1]+2, -400*x[0]]
                                                    ,[-400*x[0],200]]),
                        x0= [5,5],eps=1e-5,k=0,maxk=5000)

print(result)
